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Low-density neutron matter has been studied extensively for many decades, with a view to bet- 
ter understanding the properties of neutron-star crusts and neutron-rich nuclei. Neutron matter 
is beyond experimental control, but in the past decade it has become possible to create systems 
of fermionic ultracold atomic gases in a regime close to low-density neutron matter. In both these 
contexts pairing is significant, making simple perturbative approaches impossible to apply and neces- 
sitating ab initio microscopic simulations. Atomic experiments have also probed polarized matter. 
In this work, we study population-imbalanced neutron matter (possibly relevant to magnetars and 
to density functionals of nuclei) arriving at the lowest-energy configuration to date. For small to 
intermediate relative fractions the system turns out to be fully normal, while beyond a critical polar- 
ization we find phase coexistence between a partially polarized normal neutron gas and a balanced 
superfluid gas. As in cold atoms, a homogeneous polarized superfluid is close to stability but not 
stable with respect to phase separation. We also study the dependence of the critical polarization 
on the density. 

PACS numbers: 21.65.-f, 03.75.Ss, 05.30.Fk, 26.60.-c 



I. INTRODUCTION 

The phase structure of polarized (or population- 
imbalanced) Fermi gases is very rich. The energetically 
favored phase depends on the masses and relative popula- 
tions of the species, the interaction strength, and whether 
the interactions are long ranged. Polarized gases can ex- 
ist in a variety of physical systems, including ultracold 
trapped atomic Fermi gases [l|j, dense hadronic phases 
near the cores of neutron stars [2tJ1 , and possibly quark 
matter 

In this paper we focus on the possibility that neutrons 
in the inner crust of neutron stars are spin-polarized Q . 
Theoretical calculations suggest that this region, which is 
about a kilometer thick, features unbound neutrons while 
the protons cluster in neutron rich nuclei. While inter- 
actions with nuclei affect the properties of the unbound 
neutrons @, Q , in this work we consider an isolated neu- 
tron fluid and compute the conditions for the coexistence 
of an unpolarizcd phase and a polarized phase at densities 
roughly 30 — 250 times smaller than the nuclear satura- 
tion density. A polarized phase in the neutron star crust 
will have a larger specific heat and therefore its pres- 
ence may affect the time scales associated with observed 
crustal transient phenomena 0. Neutron superfluidity 
also affects thermal transport properties [HI hj. If the 
polarization is large enough all effects of superfluidity are 
expected to disappear. 

At the low densities we consider, the dominant term in 
the interaction between spin-up and spin-down neutrons 
is an s-wave potential with a scattering length a ~ —18.6 
fm and range r e ~ 2.7 fm. In a balanced system this at- 
traction drives the formation of Cooper pairs of spin-up 
and spin-down neutrons of opposite momenta resulting 
in a neutron superfluid whose properties can be under- 
stood qualitatively using BCS theory. All the neutrons 



arc paired and the number of up and down species is 
equal. Even at densities much smaller than the nuclear 
saturation density, a is much larger than the interparticle 
spacing and obtaining quantitatively reliable results re- 
quires non-perturbative techniques. Detailed Quantum 
Monte Carlo (QMC) simulations have calculated the en- 
ergy density £ up to densities equal to or larger than the 
nuclear saturation density fl2L [l3|. At low densities, the 
s-wave pairing gap obtained peaks at nearly 2 MeV. 

Magnetic fields tend to split the up and down Fermi 
surfaces thus disrupting the pairing and creating a polar- 
ization. For mag netic fields (B) smaller than 10 16 - 10 17 
Gauss, pairing wins (assuming no effect coming from the 
lattice of nuclei), while larger fields can create a polarized 
phase. In BCS theory, the critical field, or the resulting 
chemical potential split (<5/j, = g^^NB where 9n^n is 
the neutron magnetic moment), is called the Clogston- 
Chandrasckhar point. At this point, the BCS superfluid 
can coexist in equilibrium with a partially polarized nor- 
mal phase which is essentially a gas of free neutrons. Nu- 
merically 5/j, c = where A is the gap in the BCS 
phase. A weakly polarized gas of neutrons confined in a 
volume forms a mixed phase [l4[ consisting of the BCS 
and normal phases with volume fractions depending on 
the polarization. Hartree corrections substantially mod- 
ify the Clogston-Chandrasekhar point even for modest 
couplings Ifc^al ~ 0.5, where kp = (S^n) 1 / 3 is the Fermi 
wave vector [lfj . 

At strong coupling (mathematically given by \kpa\ > 
1) the energies of both the superfluid and the normal 
phase need to be computed non-perturbatively. In this 
paper, we perform these calculations for the normal phase 
and find the critical polarization at which it can coexist 
with the superfluid phase. We also compare the energy of 
a phase separated state with the previously evaluated @ 
homogeneous polarized superfluid phase. In this paper 
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we only consider competition between homogeneous and 
isotropic phases and ignore the possibility of promising 
but still unconfirmed scenarios (like p— wave [161 ] and in- 
homogeneous phases fl7j ) . 

As mentioned above, our microscopic ab initio calcula- 
tions may impact the observed behavior of neutron stars 
with extremely large internal magnetic fields. It is inter- 
esting to concurrently examine possible connections of 
such simulations with terrestrial experiments. Our re- 
sults may be useful for constraining Density Functional 
Theories (DFTs) that are employed to study properties 
of heavy nuclei. Equation of state results at densities 
close to the nuclear saturation density have been used for 
some time to constrain Skyrme and other other density 
functional approaches to heavy nuclei [T^]. The poten- 
tial significance of such calculations has led to a series 
of publications on the equation of state of low- de nsity 
neutron matter over the last few decades [H, fl9l - l29| . 
The density-dependence of the So gap in low-density 
neutron matter has also been used to constrain Skyrme- 
Hartree-Fock-Bogoliubov treatments in their description 
of neutron- rich nuclei [30| . Similarly, an ongoing project 
is studying the limit of extreme polarization in neutron 



matter [3 If . This limit is known in condensed-matter 
physics as the "polaron" problem: one impurity embed- 
ded in a sea of fermions. Quantities of interest here in- 
clude the polaron binding energy and its effective mass, 
which can be used to constrain Skyrme and other func- 
tionals. In this line of study, DFTs can be constrained 
by the equation of state of polarized unpaired neutron 
matter that we report here. 

Another connection of our calculations with experi- 
ment may be possible through the field of atomic ultra- 
cold Fermi gases. There, the scattering length of the in- 
teraction between two hyperfine states near a Feshbach 
resonance can be tuned by changing the magnetic field. 
Most experiments have been performed in the region 
where the scattering length is infinite (unitarity limit) 
and the effective range is negligible |32H3a | . QMC calcu- 
lations performed with the same techniques that we use 
here in the same limit give results that agree well with 



the experiments for both unpolarized [13j, l37|, |38( and po- 
larized (391 - 14^ 1 Fermi gases and this gives us confidence in 
applying the same techniques to the problem with finite 
a and r e . More directly, by tuning the scattering length 
and the effective range [43[ it may be possible to obtain a 
system with similar values of k F a and k F r e as the inner 
crust and test our results. 



II. CONSTRUCTION OF THE MIXED PHASE 

The conditions for a polarized phase (P) to be in phase 
equilibrium with an unpolarized superfluid (SF) phase 
are as follows: 

1. The pressures of the two phases should be equal. 

2. The average chemical potential of the two species 
in the phases, \x = (/i a + /Xf,)/2, should be equal. 

3. The chemical potential splitting, b\i = (p, a — ^bj/'i, 
should be less than the gap A in SF. 

We construct the mixed phase between the superfluid 
and the partially polarized normal phase (NP). The 
pressure in a polarized phase depends on the net density 
of fermions (n a + = n) and their relative fraction x = 
nb/n a where a (or f) is the majority species. 

Motivated by the unitary Fermi gas [3, |45[ , we write 
down the energy density as: 



3(6^ 2 ) 2 / 3 . . ,- 5/3 
5 2m (n ° g(a: ' n)) ■ (1) 



£(n a ,n h ) = - 
5 



For convenience we denote the energy per particle by E 
and the corresponding energy for a free gas by Efg = 
3fc|/(10m). In terms of E, g{x,n) = (E/E FG{t) ) 3 / 5 . 

For the unitary gas, the discussion of mixed phases is 
simplified greatly since the function g is only a function 
of x. In neutron matter, the interaction potential has 
inherent scales (these can be thought of simplistically as 
a and r e ). 
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To satisfy the conditions of phase equilibrium, we need to calculate the chemical potentials and the pressure, which 
include the functional dependence of g on both x and n. The chemical potential is: 

88 32/37r4/3 (iTi) 2/3 g(^") 2/3 G£r&(*>») + (!" x)^(x,n) + g(x,n)) 



on 



2 4/3 7 



(2) 



The chemical potential splitting is: 



32/3^4/3 2/3 g{Xi n) 2/3 {{ _ x _ 1} || (;Ei n) + 5(Xj n)) 



5/j, 



d£_ 

dSn 



2 4 /3 r 



(3) 



The pressure is simply P = —£(n, 5n) + /in + Sfj.Sn, where Sn = n a — rib- 

For the unpolarized superfluid (x = 1), the function g depends only on the density (5n = 0) and we have: 



(J*SF 



32/3^4/3 (zf 3 gsF ( n )W (n^(n) +g SF (n) 



2 4 /3 r 



(4) 
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To organize the discussion it is convenient to start from 
the completely polarized (x = 0) gas. Here the inter- 
actions between the only species present are weak and 
the ground state is well described by a normal phase, 
(p— wave pairing and interactions will modify the energy 
of the fully polarized neutron gas from a free Fermi gas. 
As we discuss in Sec. IIII1 we include the effect of p— wave 
interactions exactly, but ignore p— wave pairing.) 

Now consider the following thought experiment. Keep- 
ing the density n constant, change a few a fermions to b 
fcrmions so that x changes. For small x we expect that 
the phase obtained will be a partially polarized normal 
phase. As we increase the fraction of b fermions, it is 
possible that some polarized supcrfluid phase (SFP) be- 
comes the favored state of matter. Another possibility is 
that at some fractional density x c , it becomes favorable 
to form a mixed phase between the partially polarized 
normal phase (NP) and the unpolarized superfluid phase 
(SF). In the present section, we will focus on the mixed 
phase between NP and SF. In Sec. IIIII we will calculate 
the energy for the polarized normal gas (NP) phase and 
in Sec. IIVI we will look at the competition between the 
homogeneous polarized superfluid (SFP) of Ref. [f| and 
the mixed phase. 

Equating the chemical potentials and the pressures 
gives two equations in two variables, x c and Usf- n SF 
is the total supcrfluid density and tinp is the total den- 
sity of the normal phase. For x > x c one can find the 
coexistence curve in (x, n) space. Suppose that a vol- 
ume fraction v is occupied by the NP phase and the rest 
(1 — v) by the superfluid phase. Then: 



n(v) = vtlnp + (1 — v)nsF 



(5) 



x(v) 



(TT^ n NP + (l-v)^) 



(jTT c n NP + ( 1 



(6) 



where v G [0, 1]. 

The parametric dependence on v can be eliminated to 
find g(x,n,Np) along the coexistence curve. The coex- 
istence curve (for example, see dashed curve in Fig. [3]) 
should not be seen as the continuation of the NP curve 
at constant density, since the density changes along it. It 
is simply a projection of the coexistence curve in (g, x, n) 
space, on the (g, x) plane. 

To proceed with the calculation of the coexistence 
curve, we need to calculate g(x, n) for the NP phase 
and gsF(n) for the SF phase. We use Quantum Monte 
Carlo techniques to calculate the energies at a few val- 
ues of (x,n) and interpolate to determine the functional 
dependence. This is discussed next. 



III. GREEN'S FUNCTION MONTE CARLO 
SIMULATIONS 

The Hamiltonian for low-density neutron matter is: 



fe=l i<j' 



(7) 



where TV is the total number of particles. The neutron- 
neutron interaction is not purely s-wave but still some- 
what simple if one considers the AV4' formulation: [4(| 

v 4 (r) = v c (r) + v<j(r)cri ■ cr 2l (8) 

In the case of S~0 (singlet) pairs this gives: 

v s (r) = v c (r) - 3v a (r) . (9) 

However, it also implies an interaction for S=l (triplet) 
pairs: 



v P (r) = v c (r) + v a (r) 



(10) 



Ref. [471 ] explicitly included such p-wave interactions in 
the same-spin pairs (the contribution of which was small 
even at the highest density considered), and perturba- 
tively corrected the S = l,Ms = pairs to the correct 
p-wave interaction. 

In these calculations it is customary to first employ a 
standard Variational Monte Carlo simulation, which min- 
imizes the expectation value of the Hamiltonian given a 
variational wave function Vfy . At a second stage, the out- 
put of the Variational Monte Carlo calculation is used as 
input in a fixed-node Green's Function Monte Carlo sim- 
ulation, which projects out the lowest-energy eigenstate 
from the trial (variational) wave function "fy. This is 
accomplished by treating the Schrodingcr equation as a 
diffusion equation in imaginary time r and evolving the 
variational wave function up to large r. The ground state 
is evaluated from: 



*o = exp[-(i? - E T )r]* v 

~[exp[-(H - E t )At}^ x 



(11) 



evaluated as a branching random walk. The fixed-node 
approximation gives a wave function VPo that is the 
lowest-energy state with the sames nodes (surface where 
, 5 = 0) as the trial state \fv. The resulting energy Eq is 
an upper bound to the true ground-state energy. 
The ground-state energy Eq can be obtained from: 



En = 



0/ 



(*v|*o) 



QPo|ff|tto) 
(*ol*o) 



(12) 



The variational wave function is taken to be of the 
following form: 



*y(R) = n/ P (ry) [J / P (r^0n/(rv)$(R) 



(13) 
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In our earlier works we studied supcrfluid neutron matter 
and therefore used a trial wave function of the Jastrow- 
BCS form with fixed particle number. 0, [r| E3] In the 
present work we are attempting to determine the sta- 
bility of different phases. In order to do this, we have 
completely mapped out the energy of a normal (i.e. non- 
supcrfluid) neutron g function of the polarization 

at different densities. Thus, the ^(R) we used in Eq. 
(|T3| describes the particles as being in a free Fermi gas 
(i.e., lets all the correlations lie within the Jastrow func- 
tions). In other words, the wave function is composed 
of two Slater determinants (one for spin-up particles and 
one for spin-down ones): 

$ S (R) = S s (RK3s(R)jvj (14) 

where 

$s(R)at t = ^„(ri)^„(r 2 ) . . . cj>„(r Nf )] (15) 

and 

$s(R)jv; = A[<j> n {rv)<i>n{r2-) ■ ■ ■ 4>n{r N 'J\ (16) 

The primed (unprimed) indices correspond to spin-up 
(spin-down) neutrons and + N[ = N. A is the anti- 

symmetrizer and <f>n{fk) = e ik "' rfc /L 3 / 2 . 

The Jastrow part is usually taken from a lowest-order- 
constraincd-variational method calculation described by 
a Schrodinger-like equation: 

--V 2 f(r)+v(r)f(r) = A/(r) (17) 
m 

for the opposite-spin f(r) and by a corresponding equa- 
tion for the same-spin fp(r). Since the f(r) and fp(r) 
are nodeless, they do not affect the final result apart from 
reducing the statistical error. 

We calculate ground-state energies at different total 
number densities [n = (N^+N±)/L 3 ], more specifically at 
m = 6.65 x 10~ 4 , n 2 = 2.16 x 1CT 3 , and n 3 = 5.32 x 1(T 3 
fm -3 . To put these densities into perspective we can 
compare them to nuclear matter saturation density: they 
are 0.41, 1.35, and 3.32 percent, respectively, of n = 0.16 
fm~ 3 . These simulations were performed at values of 
the relative fractions chosen specifically in order to en- 
sure a full coverage of the cc-axis. More specifically, sim- 
ulations were carried out at relative fractions of x = 
0, 0.333, 0.579, 0.818, 1, corresponding to particle num- 
bers of 33 + 0, 57 + 19, 57 + 33, 33 + 27, 33 + 33, re- 
spectively. These systems are quite large (and therefore 
computationally demanding) so as to ensure a minimiza- 
tion of finite-size effects (see also Refs. (l3l l47l. |48| ) . 

The results are shown in Fig. [1] as points (circles, 
squares, and diamonds, respectively), along with cubic 
fits to the microscopic data. The latter will be useful to 
us in the following section, when we try to check the rel- 
ative stability of different phases. To facilitate the use of 
these results in connection with Eq. (jTJ) , we have divided 
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FIG. 1: (color online) Ground-state energy per particle (in 
units of the free Fermi gas energy) scaled to the power of 3/5 
for normal spin-polarized neutron matter. Shown are QMC 
results at three different total densities ni,ri2,ri3 (as points: 
black circles, red squares, blue diamonds, respectively), along 
with cubic fits to the Monte Carlo results (as lines: black solid, 
red dashed, blue dotted, respectively). Also shown using hol- 
low symbols are the results for an unpolarized superfluid. 



the ground-state energies with the energies of correspond- 
ing free spin-up Fermi gases, and raised the result to the 
3/5 power. When i = 0we see that the values are very 
close to 1, though still above it: this is due to the same- 
spin p-wave interactions. Similarly, as we increase the 
density the energy is decreased, a fact that, as we see in 
Fig. [TJ holds at any relative fraction x. The statistical 
errors of the Quantum Monte Carlo results are smaller 
than the symbols shown in the figure. Also shown to the 
right of the figure are the QMC results for an upolarized 
supcrfluid at the three densities of interest (hollow cir- 
cles, squares, and diamonds, respectively). The results 
for n\ and n-i were taken from Ref. [13]; the value at n 3 
was re-optimized for the purposes of the present work. 



IV. PHASE COMPETITION 

Having arrived at the functional dependence of the en- 
ergy of a normal neutron gas on the relative fraction 
(thereby getting information on the g(x,n) function de- 
fined in Eq. (JTJ) and appearing in the rest of Sec. [Tl]) , we 
can now use the simple cubic fits for the three densities 
m, n-i-, n-i to determine the competition between different 
phases. More specifically, Fig. [T] shows the behavior of 
the NP (normal polarized) equation of state. From the 
hollow points on the same figure we also have access to 
the energy of the SF (unpolarized superfluid) phase. Us- 
ing these two dependences, we can find out the critical 
values of the relative fraction (x c \, x c 2, Xc3, respectively) 
at the three densities we are studying. Above a critical 
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FIG. 2: (color online) The construction of the mixed phase at 
density ni . Shown is the energy of a normal neutron gas from 
Fig. [T]of the present work (black solid line), along with the 
homogeneous superfluid results of Refs. [H, [47J (red square 
points), and a tangent construction showing the phase co- 
existence curve (green dashed line). For large values of the 
relative fraction, the phase separated state is energetically 
favored with respect to a homogeneous polarized superfluid. 



relative fraction the system phase separates into a mix- 
ture of SF and NP, whereas below that value the neutron 
gas is normal, and therefore does not exhibit any of the 
widely studied features of a superfluid in the neutron star 
crust. 

To determine x c we need the functional dependence of 
g(x,n) on both x and n, which is given by the interpo- 
lations for m, ni and 713. Let us first focus on ni- In 
Eq. we also need the derivative with respect to n. We 
estimate this derivative by finite differences: (g{x, ri2) — 
g(x, ni))/(n 2 -ni) and (g(x, n 3 )-g(x, n 2 ))/ (n 3 -n 2 ). We 
estimate the derivatives for SF used in Eq.2]with a linear- 
form in \/(kpa) and kpr e . In practice, ^ is small and 
doesn't affect the value of x c and tisf significantly. We 
can then solve the conditions for coexistence and find x c 
and usf- The numerical values show slight dependence 
on the interpolation method used. A simple cubic fit as 
a function of x gives x C 2 = 0.78 and usf = 2.17 x 10~ 3 
fm -3 . A three parameter fit in (1,2;, a; 5 ' 3 ) motivated 
by 0,113 gives x c2 = 0.77 and n SF = 2.18 x 10~ 3 fm-' 3 . 
This gives some idea about the errors made in the inter- 
polation: we estimate our error bar on x c to be 0.02 for 
the three densities. This detailed procedure gives results 
very similar to a simpler procedure, namely constructing 
a tangent to the constant density NP curve for g(x,n) 
passing through the value of gsF(n) at the same density, 
which gives x C 2 = 0.78. The close agreement between 
the two approaches stems from the fact that ^ is small 
- for = 0, the coexistence curve is given exactly by 
the tangent construction (4f| . Therefore, we will use the 
tangent construction to find x c for n\ and 77,3 below. 



This analysis requires that for x < x c the normal phase 
is thermodynamically stable. That is, the eigenvalues 

,2 c 

of the matrix d " £ n are positive. For n — n 2 we can 
check and have checked explicitly that this is the case. 
Finally, we find Sii at x C 2 is given by 5[a C 2/Ef = 0.16(2) 
where Ep = (hp /2m). Comparing with AsF2/Ep ~ 
0.23 obtained from Ref. |47|, we see that A$F2 > <^Mc2 
satisfying the third condition for coexistence. 

In Fig. [3] we have plotted as a black solid line the 
n\ results from Fig. [TJ at this density the perturbative 
5 = 1, Mg = correction in the Hamiltonian men- 
tioned in Sec. IIIII is the smallest, leading to the high- 
est degree of confidence in the accuracy of our micro- 
scopic results. We focus on large values on the a;-axis for 
reasons that will soon become clear. We also show as 
red square points the unpolarized superfluid result from 
Ref. [43] (at x = 1) and the homogeneous polarized su- 
perfluid results of Ref. @ (at smaller x). In addition, 
we show a coexistence curve which is arrived at by a 
tangent construction from the unpolarized superfluid to 
the norm al p olarized curve when the y-axis is chosen as 
we have [45||. This coexistence line meets the NP curve 
at x c \ = 0.74 (Sfid/Ep = 0.2). As can be clearly seen 
in this Figure, the coexistence line lies below the homo- 
geneous polarized superfluid Quantum Monte Carlo re- 
sults, implying that these are not energetically favored 
(foreseeing this possibility, Ref. @ included a critical 
assumption explicitly excluding the possibility of phase 
separation). This situation follows the behavior of ultra- 
cold fermionic gases at unitarity: [H,|!qI,IH! there the crit- 
ical value of the relative fraction was closer to x c = 0.44. 
The value we find is different for a variety of reasons: 
first, ni corresponds to kpa = —5, which is somewhat 
smaller than kpa = 00. More importantly, the neutron- 
neutron interaction is characterized by a finite effective 
range and the presence of higher partial waves. The lat- 
ter point explains why when we repeat this exercise for 
71.2 {kpa = —7.4) and 77.3 (kpa = —10) the correspond- 
ing results (from tangent constructions) for the critical 
fraction are x C 2 = 0.78 (Sfi^/Ep = 0.16) and x c s = 0.88 
{&Hcz/Ep = 0.07), respectively. The trend appears to 
be that at larger density (stronger coupling) the critical 
fraction is moving toward 1, implying that farther down 
toward the core of the neutron star polarized neutron 
matter quickly becomes normal for the vast majority of 
possible polarizations. 

The value of 5ii c is relevant for neutron star phe- 
nomenology since whether 5fj, is greater or less than 8fi c 
determines the phase. This is seen more simply in the 
grand canonical ensemble where we can compare the 
pressure as a function of the chemical potentials. In Fig. [3] 
we plot the pressure as a function of S[i for fixed /i = 1.94 
MeV. For this value of 11 the density of the superfluid 
phase at the critical 5/i is 712. For <5/i < <5/i c the pressure 
of the SF phase is larger and the system is unpolarised. 
Even if the magnetic field is not large enough to polarize 
the system, certain thermal obscrvablcs, for example the 
specific heat, are sensitive to the difference between 8ti 
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FIG. 3: (color online) Pressure as a function of S/j, for fixed 
fi = 1.94 MeV. For Sfi > 5fi c = 0.53 MeV the unpaired phase 
(bold black curve) has a higher pressure. At S[i c [marked by 
the vertical line (brown online)], the SF and the NP phases 
can coexist. For smaller Sfj,, the pressure of the SF phase 
marked by the horizontal line (blue online) is larger than the 
pressure of the NP phase marked by the dashed line (red 
online), and the SF phase is favored. On the other extreme 
for Sfi = 4.18 MeV, x = 0.0005. We choose fi so that the 
density of the normal phase at 8fi c is rii- 



and 8fi c . 

Determining the critical relative fractions is also rel- 
evant: if a magnetic field is strong enough to polarize 
neutron matter sufficiently (x < x c ), the gas ceases to be 
superfluid even at zero temperature. At such large mag- 
netic fields, processes that rely on the presence of a bulk 
superfluid in the inner core of neutron stars, for example 
a recently proposed heat-conduction mechanism which 
requires superfluid phonons [lj|, would no longer be op- 
erational. On a different note, a precise determination 
of these fractions (along with the rest of the equation 
of state) is also significant for Skyrme functional prac- 
titioners: for the region where the NP phase is the true 
ground state of the system, the results shown in Fig.[T]are 
a microscopic constraint to nuclear energy-density func- 
tional theory, which can help improve its predictions on 
neutron-rich heavy nuclei. 



V. SUMMARY & CONCLUSIONS 

In summary, we have studied spin-polarized low- 
density neutron matter at many values of the polariza- 
tion. Our Quantum Monte Carlo approach provides tight 
upper bounds which are expected to be quite close to 
the true ground state-energy of the system. At very-low 



density the Hamiltonian becomes quite simple and we 
include its dominant well-known terms. We have calcu- 
lated the energy as a function of polarization using the 
AV4' interaction at three different densities. 

One of our main results is a determination of the criti- 
cal relative fraction values above which the system phase 
separates into an unpolarized superfluid and a polarized 
normal neutron gas. Above these values, a homogeneous 
polarized superfluid is found not to be energetically fa- 
vored. Below these values, the neutron gas is completely 
normal. Both these facts can lead to astronomically im- 
portant consequences. Importantly, the trend of our re- 
sults seems to imply that at slightly larger densities x c 
would be very close to one: this would mean that even at 
zero temperature a small polarization would be enough 
to close the superfluid pairing gap, and thus easily lead 
to larger polarization. In a non-astrophysical context, it 
is conceivable that our results could be tested directly by 
using ultracold fcrmionic atom gases with unequal spin 
populations and a finite effective range. In cold atoms, 
Quantum Monte Carlo simulations of spin-polarized mat- 
ter have in recent times been repeatedly used as input to 
computationally less demanding density-functional the- 
ory approaches. Similarly, we expect that the results 
presented in this work can also be used as input to self- 
consistent mean- field models of nuclei. 

A possible direction for future work, that would be 
interesting to both nuclear astrophysicists and Skyrme 
practitioners, would be a study of the behavior of po- 
larized and unpolarized superfluid and normal neutron 
matter in an external periodic potential. Such a study of 
the static response of neutron matter could also in prin- 
ciple be guided (or guide) analogous research related to 
optical lattice experiments with cold atoms. 

In the neutron star context, it will be useful to analyze 
the finite temperature phase diagram to understand how 
the specific heat is modified as a function of the chemi- 
cal potential splitting, and how this modification affects 
transient phenomena in the neutron star crust. 
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